Número de SADs preditas não refutadas
- Padrões gerais
- GAMM binomial: número de predições não refutadas
diff_S = (S_obs - S_MN)/S_obs
- Padrão Geral
- GAMM: diff_S
  - Modelo Cheio e Gráficos Diagnosticos
    - Observado e Predito
U - taxa de especiação necessária para obter a riqueza observada no equilíbrio
- Padrões Gerais
- GAMM
Parâmetros dos Modelos Neutros

Número de SADs preditas não refutadas

Padrões gerais

Figura 1.1 logito(número de SADs não refutadas) ~ p * k * MN. A linha azul é uma estimativa baseada em ‘loess’.

A figura 1.1 mostra que o padrão geral dos dados não é linear na escala da função de ligação, assim utilizamos GAM. Possibilitamos até 1 intercepto por modelo neutro por sítio de amostragem (MN|SiteCode); para modelos com variáveis de dispersão continua não foi possível um smoother por Sítio de Amostragem e Modelo Neutro

GAMM binomial: número de predições não refutadas

Comparação de Modelos Cheios

l_md <- vector("list",6)
names(l_md) <- c("k 1|SiteCode", "k MN|SiteCode",
                 "k_1 1|SiteCode","k_1 MN|SiteCode",
                 "d/L 1|SiteCode","d/L MN|SiteCode")
# muitos parametros para poucos dados "k_1 k_1 * MN|SiteCode","d_Lplot d_Lplot * MN|SiteCode")
# k
l_md[[1]] <- gam(cbind(n_nRef,100-n_nRef) ~ k * MN + s(p.z,k,by=MN,bs="fs") +
                   s(SiteCode,bs="re"),
                   family = "binomial",
                   data = df_resultados, method = "REML")
l_md[[2]] <- gam(cbind(n_nRef,100-n_nRef) ~ k * MN + s(p.z,k,by=MN,bs="fs") +
                   s(SiteCode,by=MN,bs="re"),
                   family = "binomial",
                   data = df_resultados, method = "REML")
# k_1
l_md[[3]] <- gam(cbind(n_nRef,100-n_nRef) ~ MN + s(p.z,by=MN,bs="tp") + s(k_1.z,by=MN,bs="tp") + ti(p.z,k_1.z,by=MN,bs=c("tp","tp")) +
                   s(SiteCode,bs="re"),
                   family = "binomial",
                   data = df_resultados, method = "REML")
l_md[[4]] <- gam(cbind(n_nRef,100-n_nRef) ~ MN +  s(p.z,by=MN,bs="tp") + s(k_1.z,by=MN,bs="tp") + ti(p.z,k_1.z,by=MN,bs=c("tp","tp")) +
                   s(SiteCode,by=MN,bs="re"),
                   family = "binomial",
                   data = df_resultados, method = "REML")
# d/L
l_md[[5]] <- gam(cbind(n_nRef,100-n_nRef) ~ MN + s(p.z,by=MN,bs="tp") + s(d_Lplot.z,by=MN,bs="tp") + ti(p.z,d_Lplot.z,by=MN,bs=c("tp","tp")) +
                   s(SiteCode,bs="re"),
                   family = "binomial",
                   data = df_resultados, method = "REML")
l_md[[6]] <- gam(cbind(n_nRef,100-n_nRef) ~ MN +  s(p.z,by=MN,bs="tp") + s(d_Lplot.z,by=MN,bs="tp") + ti(p.z,d_Lplot.z,by=MN,bs=c("tp","tp")) +
                   s(SiteCode,by=MN,bs="re"),
                   family = "binomial",
                   data = df_resultados, method = "REML")
AICctab(l_md,weights=TRUE)

##                 dAICc   df               weight
## k MN|SiteCode       0.0 567.902990229433 1     
## k_1 MN|SiteCode 10456.4 255.87941997766  <0.001
## d/L MN|SiteCode 14414.2 255.878481862603 <0.001
## k 1|SiteCode    29130.6 480.710534787585 <0.001
## k_1 1|SiteCode  38331.6 162.939249908097 <0.001
## d/L 1|SiteCode  43420.3 162.955609731983 <0.001

# save(l_md,file="~/Documentos/Doutorado/artigo_mestrado/1A.P_MOVER/Support Information/md_nRef.Rdata")
# load("~/Documentos/Doutorado/artigo_mestrado/1A.P_MOVER/Support Information/md_nRef.Rdata")

O único modelo plausível considera a variável k (fator) e estrutura aleatória MN|SiteCode.

##       dAICc df               weight
## by=MN   0.0 567.902990229433 1     
## by=k   47.0 497.198353791482 <0.001

Diagnostico do Modelo Cheio

Figura 1.2 appraise(gam( cbind(n_nRef,100-n_nRef) ~ k*MN + s(p.z,MN,by=k), binomial(logit) )

O modelo mais plausível não apresenta boa congruência com os dados: a distribuição binomial não parece se adequar bem aos dados e o predito e observado podem ser muito divergentes.

Figura 1.3 Probabilidade de não refutar uma SAD neutra para cada modelo estatístico (título dos paineis); na esquerda o modelo mais plausível, na esquerda o segundo modelo plausível (>1000 deltaAIC). Os valores médios de cada modelo neutro apresentam semelhante entre os dois gráficos, contudo o modelo mais plausível apresenta maior flexibilidade na descrição dos modelos.

Figura 1.4 Observado e predito pelo modelo mais plausível. Eixo x = prop. de cobertura vegetal na paisagem; Eixo y: número de SADs não refutadas. Em branco o padrão geral estimado e em vermelho o erro padrão.

o modelo não apresenta bom ajuste aos dados:
i) MNEE: de maneira geral a tendência global parece estar de acordo com o observado. Contudo, para k<0.5, a estimativa parece ser muito influenciada por pontos extremos.
ii) MNEI: para todo k, a tendência global não parece apresentar boa congruência com o observado principalmente quanto aos valores médios observados, e.g., para k>0.50 a tendência global parece subestimar o observado. Para k<0.5 há um ponto que parece ser de grande influência (k=0.25 p->0.5).

Figura 1.5 residuals.gam(l_md[[“k MN|SiteCode”]]) ~ p (~k*MN).

Os resíduos do tipo ‘deviance’ não indicam erro na estimativa quanto à tendência global; o quê para mim está errado principalmente para MNEI.

Melhoria do modelo mais plausível: CONSTRUÇÃO

Opções

Comparação metodo para estimar sp e classe de bs

smoother type

##    dAICc  df               weight
## cr    0.0 576.311752135652 1     
## cc 2939.1 543.980281194929 <0.001
## tp 3290.1 567.902990229433 <0.001
## ps 3506.6 559.143708020173 <0.001
## cp 6113.9 566.50730949499  <0.001

##         dAICc df               weight
## GACV.Cp   0.0 582.038923755289 0.5   
## GCV.Cp    0.0 582.038923755289 0.5   
## REML     38.0 576.311752135652 <0.001
## ML       38.4 576.287674101362 <0.001

Figura 1.6 Predito pelo gam com method=GACV.Cp e bs=cr. O intervalo de confiança foi removido por apresentar valores negativos.

Apesar das mudanças na configuração o modelo continua estimando MNEI errado.

Resumo dos achados da comparação:

GCV.Cp e GACV.Cp são os métodos mais plausíveis para o smoothing type padrão (“tp”) e para o mais plausível (“cr”).
Porém estes métodos de estimativa do parâmetro de smoothing apresentam ajustes estranhos aos padrões observados.
O smoother type “mais plausível”cr" ainda apresenta o mesmo problema do anterior: MNEI não está fazendo bom ajuste.

Como não consegui configurar os métodos GCV.Cp e GACV.Cp para um padrão mais adequado sigo com o terceiro mais plausível (“REML”, o padrão).
Vou explorar ajustar modelos separados para cada conjunto de dados

funções de ligação para distribuição binomial

##         dAICc  df               weight
## cloglog    0.0 572.962297266108 1     
## logit   2249.8 576.311752135652 <0.001
## probit  3640.8 575.395171500865 <0.001

Figura 1.7 Gráficos Diagnósticos do modelo mais plausível binomial(link=cloglog)

A qualidade do modelo não melhorou.

Figura 1.8 Predição pelo modelo com função de ligação “cloglog”

Para MNEI o problema persiste.

transformação do conjunto de dados e distribuicao normal e gamma

##               dAICc  df               weight
## normal           0.0 347.433185775591 1     
## gamma log     4941.1 288.45300495345  <0.001
## gamma inverse 5271.7 264.637855151765 <0.001

Figura 1.9 Resíduos logit_nRef normal

Figura 1.10 Predito pelo modelo com distribuição normal

Apesar de possibilitar uma boa estimativa da média o modelos é muito

Avaliação do conjunto de dados separadamente

MNEE

##                        dAICc  df               weight
## k 1|SiteCode logit        0.0 261.404710989706 1     
## k 1|SiteCode probit     663.6 261.341234025727 <0.001
## k 1|SiteCode cloglog   1254.7 256.432441274999 <0.001
## d/L 1|SiteCode logit   7212.7 127.452558167134 <0.001
## d/L 1|SiteCode probit  7212.7 127.452558167134 <0.001
## d/L 1|SiteCode cloglog 7212.7 127.452558167134 <0.001
## k_1|SiteCode cloglog   7865.4 127.890778649655 <0.001
## k_1|SiteCode probit    8539.6 127.940128196497 <0.001
## k_1|SiteCode logit     8650.6 127.878762400437 <0.001

Figura 1.11 Gráficos Diagnósticos do modelo mais plausível para o conjunto de dados MN:EE

Figura 1.12 Predito para o conjunto de dados MNEE

Ainda há overfitting para alguns pontos.

MNEI

##                        dAICc   df               weight
## k 1|SiteCode cloglog       0.0 255.938425616834 1     
## k 1|SiteCode logit      3522.0 265.551522242843 <0.001
## k 1|SiteCode probit     4252.3 272.189573651799 <0.001
## k_1|SiteCode cloglog    6003.1 127.972940662994 <0.001
## k_1|SiteCode logit     10150.9 127.977463144225 <0.001
## k_1|SiteCode probit    10727.2 127.984432980891 <0.001
## d/L 1|SiteCode logit   15681.2 127.9410891041   <0.001
## d/L 1|SiteCode probit  15681.2 127.9410891041   <0.001
## d/L 1|SiteCode cloglog 15681.2 127.9410891041   <0.001

Figura 1.13 Gráficos Diagnósticos do modelo mais plausível para o conjunto de dados MN:EI

Figura 1.4 Predito para o conjunto de dados MNEI

diff_S = (S_obs - S_MN)/S_obs

Padrão Geral

Figura 2.1 diff_S = (S_obs - S_MN)/S_obs + (-1) * min(diff_S) + 0.01, inclinação_MNEE ~ 0 para todo k.

Parece que existem alguns outliers. De maneira geral MNEE apresenta boa congruência com os dados enquanto o comportamento de MNEI depende da dispersão: para k->100% o desvio aumento para ps extremos (bimodal/quadratica), apresentando boa congruência para proporção intermediárias de habitat, mas sempre subestimando S_obs. Com o aumento da capacidade de dispersão as estimavas de S se tornam mais próximas ao observado para p>0.5, superestimando o observado para valores elevados de p e k>0.5 o que leva a reduzir a porção de p em que MNEI faz uma boa aproximação.

GAMM: diff_S

Para modelar o padrão não linear dos dados vou utilizar GAMM: normal, Gamma(log) e Gamma(inverse) com duas estruturas aleatórias (1|SiteCode e MN|SiteCode).

Modelo Cheio e Gráficos Diagnosticos

##                           dAICc  df               weight
## normal id MN|SiteCode        0.0 204.716290217321 1     
## normal id 1|SiteCode      3500.4 201.919668132722 <0.001
## Gamma inverse MN|SiteCode 7545.2 183.096851367051 <0.001
## Gamma log MN|SiteCode     7899.2 183.888877817823 <0.001
## Gamma inverse 1|SiteCode  8599.7 172.829381805835 <0.001
## Gamma log 1|SiteCode      8730.6 169.279594585853 <0.001

Figura 2.2 Graficos diagnosticos do modelo mais plausível para diff_S.

Parece que o pressuposto de distribuição normal não é adequado, assim como parece haver alguma tendência nos dados que não está sendo detectada. Especulo que existam outliers que estão inflenciando muito a estimativa dos dados.

Observado e Predito

Figura 2.3 Diferença na riqueza observada e estimada pela proporção de cobertura vegetal.

Os modelos ajustados para MNEI estão sistematicamente com problemas de errar a média (figura 1.1 e 1.4; figura 2.1 e 2.3). Especulo que talvez alguns pontos estão influenciando muito as estimativas centrais ou estou usando configuração errada nas funções. Além disso para MNEI percebo que a tendência geral não está descrevendo o padrão dos dados, em especial para k->1.

U - taxa de especiação necessária para obter a riqueza observada no equilíbrio

Padrões Gerais

Figura 3.1 Padrões gerais de U_med: ~ k (% de propágulos até a vizinhança imediata da planta progenitora); ~ Stotal (riqueza observada na área amostral). E o padrão empírico encontrado nos dados Stotal ~ p

GAMM

Como os modelos com estrutura aleatória mais complexa levam muito tempo para estimar, irei começar com a estrutura aleatória mais simples para avaliar qual a melhor variável de dispersão e se for necessário refazer a seleção da estrutura aleatória.

Tabela de Seleção do Modelo Mais plausível:

##                    dAICc  df               weight
## d log gamma           0.0 121.264568083601 1     
## 1-k log gamma       174.1 125.474107334272 <0.001
## k.f log gamma       227.4 154.787008201268 <0.001
## d log normal        394.5 124.471358521882 <0.001
## d inverse gamma     766.2 127.68541117751  <0.001
## 1-k log normal      855.7 127.94024990325  <0.001
## k.f log normal      889.4 170.667471702588 <0.001
## 1-k inverse gamma   983.2 124.456458934862 <0.001
## k.f inverse gamma  1002.8 140.80259373167  <0.001
## d inverse normal   1136.4 121.957626137592 <0.001
## k.f inverse normal 1345.5 165.620770162204 <0.001
## 1-k inverse normal 1354.2 122.442084920132 <0.001
## 1-k id normal      2163.0 125.186614559034 <0.001
## d id normal        2189.1 123.301338427342 <0.001
## k.f id normal      2242.8 160.04673943848  <0.001

Segue gráficos diagnosticos do modelo mais plausível

Figura 3.2 Gráficos Diagnóstico do modelo mais plausível

O modelo apresenta boa congruência com os dados, porém parece existir outliers no entando.

Figura 3.3 Métricas de influência dos dados (broom::augment())

Exclui os dois sítios com valores de grande influência do conjunto de dados:

Figura 3.4 Grafico diagnostico do modelo mais plausivel sem os outliers.

Figura 4.5 Efeitos Parciais do modelo mais plausível

Parâmetros dos Modelos Neutros

d / L_plot

padrões gerais

Figura 5.1.1 d/Lplot ~ k

Descrição Estatística

## [1] "exp(coef):"

## (Intercept)       k0.95        k0.9       k0.85        k0.8       k0.75 
##   0.6715478   1.4153327   1.7360273   2.0081002   2.2669913   2.5268653 
##        k0.7       k0.65        k0.6       k0.55        k0.5       k0.45 
##   2.7936350   3.0726860   3.3755022   3.7066527   4.0756037   4.4953275 
##        k0.4       k0.35        k0.3       k0.25        k0.2       k0.15 
##   4.9791588   5.5515868   6.2545191   7.1436593   8.3360091  10.0630874 
##        k0.1       k0.05 
##  12.9390064  19.3670608

Figura 5.1.2 Gráficos Diagnósticos do modelo gam(d ~ k + s(DA,by=k), family=Gamma(link=“log”)). Exponenciação dos coeficientes parametricos estimados pelo modelo.

m’

padrões gerais

Figura 5.2.1 m’ ~ p (~k). A variável m’ é uma função de p, d e I(1/J)

Descrição Estatística

##        dAICc  df               weight
## model2    0.0 158.703955395811 1     
## model3  488.0 43               <0.001
## model1 7043.7 113.502207953068 <0.001

Figura 5.2.2 Diagnostico do modelo mais plausivel

theta

padrões gerais

Figura 5.3.1 Padrão Geral de theta considerando escala log

Descrição Estatística

–>

–> –> –> –>

–>

–> –> –> –>

–>

–> –> –> –> –>

–> –> –> –> –> –>

–> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –>

–>

–> –> –> –>

–>

–> –> –> –>

–> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –>

–>

–> –> –> –>

–>

–> –> –> –> –> –> –> –>

–> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –> –>

–>

Estatistica

Mori, Danilo

31/10/2019

Número de SADs preditas não refutadas

Padrões gerais

GAMM binomial: número de predições não refutadas

Comparação de Modelos Cheios

Diagnostico do Modelo Cheio

Melhoria do modelo mais plausível: CONSTRUÇÃO

Comparação metodo para estimar sp e classe de bs

funções de ligação para distribuição binomial

transformação do conjunto de dados e distribuicao normal e gamma

Avaliação do conjunto de dados separadamente

MNEE

MNEI

diff_S = (S_obs - S_MN)/S_obs

Padrão Geral

GAMM: diff_S

Modelo Cheio e Gráficos Diagnosticos

Observado e Predito

U - taxa de especiação necessária para obter a riqueza observada no equilíbrio

Padrões Gerais

GAMM

Parâmetros dos Modelos Neutros

d / L_plot

padrões gerais

Descrição Estatística

m’

padrões gerais

Descrição Estatística

theta

padrões gerais

Descrição Estatística